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ABSTRACT 


A technique for approximating the modal aerodynamic influence coefficient (AIC) matrices 
by using basis functions has been developed. A process for using the resulting approximated 
modal AIC matrix in aeroelastic analysis has also been developed. The method requires the 
unsteady aerodynamics in frequency domain, and this methodology can be applied to the 
unsteady subsonic, transonic, and supersonic aerodynamics. The flutter solution can be found by 
the classic methods, such as rational function approximation, k, p-k, p, root locus et cetera. The 
unsteady aeroelastic analysis using unsteady subsonic aerodynamic approximation is 
demonstrated herein. The technique presented is shown to offer consistent flutter speed prediction 
on an aerostructures test wing (ATW) 2 and a hybrid wing body (HWB) type of vehicle 
configuration with negligible loss in precision. This method computes AICs that are functions of 
the changing parameters being studied and are generated within minutes of CPU time instead of 
hours. These results may have practical application in parametric flutter analyses as well as more 
efficient multidisciplinary design and optimization studies. 
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INTRODUCTION 


Supporting the Aeronautics Research Mission Directorate guidelines, the National 
Aeronautics and Space Administration (NASA) Dryden Flight Research Center (DFRC) 
(Edwards, California) is developing an object-oriented multidisciplinary design, analysis, and 
optimization (MDAO) tool (ref. 1). This tool will leverage existing tools and practices, and allow 
the easy integration and adoption of new state-of-the-art software. At the heart of the object- 
oriented MDAO tool is the central executive module (CEM) as shown in figure 1 . In this module, 
the user will choose an optimization methodology, and provide starting and side constraints for 
continuous as well as discrete design variables and external file names for interface variables, 
which communicate between the CEM and each analysis module. The structural analyses 
modules such as computations of the structural weight, stress, deflection, buckling, and flutter 
and divergence speeds have been developed and incorporated into the object-oriented MDAO 
framework. 

In general, obtaining aerodynamic influence coefficients (AIC) by direct calculation for 
integration into preliminary design activities involving disciplines such as aeroelasticity, 
aeroservoelasticity, and optimization is, at present, a costly and impractical venture (ref. 1). With 
the increasing complexity of the configuration and the increasing fidelity of the aerodynamic 
equations, the computational costs increase rapidly. The unsteady aeroelastic analysis and design 
optimizations are a challenging task. The time required for unsteady computations whether in 
frequency domain or time domain will considerably slow down the whole design process. Also, 
these analyses are usually performed repeatedly to optimize the final design. Even though the 
computational costs may be reduced by the use of advanced algorithms and improved computer 
hardware processing speeds, these Ml aeroelasticity analyses cannot be incorporated effectively 
within a preliminary design and optimization environment. For example, using the ZAERO™ 
(ZONA Technology Incorporated, Scottsdale, Arizona) code to generate the modal AIC matrix 
for the IKHANA (a modified Predator® B, General Atomics Aeronautical Systems, Incorporated, 
San Diego, California) aircraft with existing AIC data, one Mach number, and 16 reduced 
frequencies, takes an average of 30 min on an Intel® Core™2 2.80 GHz CPU computer (Intel 
Corporation, Santa Clara, California) (ref. 1). One of the examples presented in this paper, the 
hybrid wing body (HWB) aircraft, takes about 50 min to generate the modal AICs with existing 
AIC data. 

As a result, there is considerable motivation to be able to perform aeroelastic calculations 
more quickly and inexpensively using basis function approximation (BFA) method. One of the 
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primary goals behind the current development is to reduce the computation time for generating 
modal AIC matrices during the optimization procedure. 

BASIS FUNCTION APPROXIMATION 

The BFA method requires the unsteady aerodynamics to be represented in the frequency 
domain. In this study, the BFA for the subsonic and supersonic speeds are discussed. The flutter 
solution can be found by the classic methods, such as rational function approximation, k, p-k, p, 
root locus et cetera. A process that efficiently incorporates an approximated modal AIC matrix 
into the MDAO tool for design optimization at a reasonable computational cost has been 
developed and is outlined in the flowchart shown in figure 2. 

In linear algebra, a basis is a set of vectors that in a linear combination can represent every 
vector in a given vector space or free module, and no element of the set can be represented as a 
linear combination of the others (ref. 2). In other words, a basis is a linearly independent spanning 
set. 


This paper discusses an effective approach for approximating a modal AIC matrix for 
optimization with flutter speed constraints. Consider the steps of the approximation process 
depicted in the flowchart given in figure 2. Steps 1 and 2 are performed only once before starting 
the optimization process. Steps 3, 4, and 5 are done iteratively for aeroelasticity calculation. 

In step 1, a set of representative basis functions T is defined and intended to capture salient 
features of the modal responses the airplane is expected to encounter in the design space. 
Structural mode shapes obtained from different mass or various stiffness configurations of the 
airplane can be used as the basis functions. These basis functions are comparison functions 
(ref. 3) since all the geometric and natural boundary conditions of the airplane are satisfied. Each 
mode shape of the airplane with the target configuration can be approximated as a linear 
combination of a set of the basis functions. The target configuration in this paper is simply an 
arbitrary design point within the design space and is used as a check case. 

In step 2, the representative basis AIC matrices Q are computed corresponding to the 
representative basis functions defined in step 1 at any Mach number and reduced frequency. 
These basis AIC matrices are used as input for approximate modal AIC matrix calculation in step 
4. 


In step 3, for a set of given design structural mode shapes 0 , each mode shape is decomposed 
in a linear combination of the basis functions. The z'-th target mode shape vector 0 (z'=l, 2, ..., m, 

where m is the number of mode shapes) is approximated through the use of a least squares 
method of the representative basis functions as shown in equation (1), 


n 

0 / ® X PkVk 


k = 1 


( 1 ) 
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where a vector i//£ is the k - th basis function and a scalar coefficient /3 l , is the modal participation 
factor of the k - th basis function on the /-th mode shape; n is the number of basis functions. 

In step 4, an approximate modal AIC matrix Q is computed based on a basis AIC matrix Q 
and modal participation factors. Let matrices Q, A, and O be a modal AIC matrix, a general AIC 
matrix, and a modal matrix, respectively. Then the modal AIC matrix Q can be defined as shown 
in equation (2), 


Q = <fr r A<& 


where O = [0, d 2 ■■■<!>,„ \ 


( 2 ) 


Therefore the /-th row and /-th column element of the modal AIC matrix Q can be written as 
shown in equation (3), 


Qij = 4>I Mj 


( 3 ) 


Substituting equation (1) into equation (3) is given in equation (4), 


Q ii= 


ij~Yi ™rj 


V fc=l 


Z Pivk 


= Qn 


\k = 1 


Qij = Z Z P'sPlwl A Vr 


r=\ 5=1 


(4) 


The detailed derivations of equation (4) are listed in Appendix A. Let the basis AIC matrix Q be 
defined as shown in equation (5), 


Q = '¥ t A'¥ 


( 5 ) 
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then the 5-th row and r-th column element of the Q can be written as shown in equation (6), 


Qsi — Vs A-Vr 


( 6 ) 


Substituting equation (6) into equation (4) gives equation (7) for the approximate modal AIC 
matrix. 



iimQsr 


r=ls=l 


(7) 


The basis AIC matrix Q in equation (5) can be computed and saved before starting 
optimization as mentioned in step 2. During optimization, the mode shapes, b/ (7=1 , 2, ..., m), are 

fitted using basis functions i//^ and modal participation factors /3[_ as shown in equation (1), and 
the approximate modal AIC matrix Q can be computed using equation (7). This approximate 
modal AIC matrix will be used for flutter analysis in step 5. The resulting flutter speed and 
frequency could be used for optimization. In step 5, the flutter solution can be found by the 
classic methods, such as rational function approximation, k, p-k, p, root locus et cetera. 

APPLICATIONS 

In order to validate the proposed BFA technique in subsonic flight regimes, the flutter results 
using this approximate method are compared with the direct flutter results. First a simple model 
from the aerostructures test wing (ATW) 2 program was chosen. Subsequently, a complex 
practical problem HWB model was analyzed. In this study, structural mode shapes obtained from 
different mass configurations of the airplane are used as basis functions. 

A Modified Aerostructures Test Wing 2 

The proposed technique has been applied and validated using a modified ATW2. The original 
ATW2 (ref. 4) test article was actually designed, built, and flight-tested at NASA DFRC as 
shown in figure 3. This wing was cantilevered from a center station pylon and flown on the 
McDonnell Douglas (now The Boeing Company, Chicago, Illinois) NF-15B test bed aircraft. 

In this study, calculated structural mode shapes obtained from different mass configurations 
of the wing are used as basis functions. To create each different mass configuration, the mass at 
each mass point is represented using an additional concentrated mass element attached through 
rigid body elements to the finite element model as shown in figure 4. These alterations 
unintentionally increased the effective wing stiffness and the first natural frequency over that of 
the actual ATW2 configuration (ref. 5) in spite of the mass increase because of the rigid elements. 

The total weight and center of gravity of the various fictitious and target ATW2 mass 
configurations are listed in table 1. The direct flutter analysis of the target configuration (mass = 
0.0 units) is performed at Mach 0.82 with 16 reduced frequencies using MSC/NASTRAN (MSC 
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Software Corporation, Santa Ana, California) and ZAERO™ codes (refs. 6-7). In order to have 
better flutter solution convergence for flutter speed, the first 10 structural modes are included in 
the calculation. While preparing for the approximate flutter solution, three different mass 
configurations, 0.3, 0.6, and 0.9 units were selected for generating the basis functions. 

In order to be able to capture more of the target mode shapes, larger basis functions were 
generated. Fifteen structural mode shapes for each mass configuration were generated using 
MSC/NASTRAN. The natural frequencies of the various fictitious and target ATW2 mass 
configurations are listed in table 2. According to Pak and Lung (ref. 5), the primary natural modes 
for the first flutter mode are mode numbers 1, 2, and 3; and all of the higher modes are the 
secondary modes. The first three natural frequencies and mode shapes are the first bending, the 
first torsion, and the second bending modes, respectively and are shown in figure 5. A total of 45 
basis functions were generated. Then the 45 by 45 basis AIC matrix corresponding to those 45 
basis functions was computed for each of 16 reduced frequencies at Mach 0.82 using ZAERO™. 
The mode shapes of the target configuration are fitted using the basis functions together with a 
least squares method; and then the approximate modal AIC matrices were computed, and the 
flutter analysis was performed. 

For the present test case, ATW2 mode shapes, generated from direct method and BFA 
method, are plotted in figures 5a, 5b, and 5c. The approximated mode shapes are very well 
matched to the target mode shapes. The difference of the first three mode shapes between these 
two methods is plotted in figure 6. 

The results of the matched flutter analysis using the direct method and the BFA method are 
summarized in table 3. Three percent damping was used for the flutter speed and frequency 
computation. Both methods predicted the same flutter speed of 543.59 KEAS at a flutter 
frequency of 49.01 Hz. The speed versus damping, V-g, and speed versus frequency, V-co curves 
from the direct and BFA methods at Mach 0.82 are given in figures 7 and 8. 

Table 4 lists the computational cost of computing the modal AIC with a given structural 
mode shape at one Mach number using the direct method and the BFA method. The comparison 
was done on an Intel® Core™2 2.80 GHz CPU computer. The computational cost of generating 
the basis functions before the parametric study is not included in the comparison. 

Hybrid Wing Body Vehicle 

The modified ATW2 is a simple cantilevered wing model and a straight forward application. 
A hybrid wing body (HWB) aircraft (ref. 8) was selected for a second and more challenging 
demonstration of the BFA method. 

In this HWB application, the basis functions are created based on three different total fuel 
conditions of the vehicle, configuration number 1, 2, and 3. The target model weight is bounded 
by the range of the weight of configurations 2 and 3. In the finite element model, the various fuel 
weights are modeled using a concentrated mass element attached directly to each of the node 
points in the fuel tank area as shown in figure 9 (ref. 9). The total weight and center of gravity 
location of the different fuel weight and target configurations are listed in table 5. 

With target configuration, the direct flutter analysis is performed at Mach 0.50, with 40 
structural modes and 14 reduced frequencies using MSC/NASTRAN and ZAERO™ codes. A 
classical wing bending torsion type of flutter was predicted in this HWB example. The first 
fourteen flexible natural frequencies and the first four flexible mode shapes are listed in table 6 
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and figure 10, respectively. The first four modes are the first symmetric and anti-symmetric 
bending, and the first symmetric and anti-symmetric torsion. 


For the approximate flutter solution, a total of 150 basis functions are generated from 3 
different fuel configurations and 50 structural mode shapes for each fuel configuration using 
MSC/NASTRAN. Then the 150 by 150 basis AIC matrix corresponding to those 150 basis 
functions is computed for each of 14 reduced frequencies at Mach 0.50 using ZAERO™. 

The flutter boundary of the HWB at Mach 0.50 using BFA method is compared to that from 
the direct method. The flutter boundaries are summarized in table 7. The flutter speed is 
normalized with respect to direct method flutter speed. With 3-percent structural damping, the 
normalized flutter speed using the BFA method is 0.991, and the flutter frequency is 2.671 Hz. 
The percentage errors of the flutter speed and frequency are 0.86 percent and 0.15 percent. The 
speed versus damping, V-g, and the speed versus frequency, V-co, curves of the HWB at Mach 
0.50 obtained from the direct and BFA methods are given in figures. 11 and 12. The BFA offers 
consistent flutter speed and frequency predictions on the HWB target configuration. See table 4 
for the computational cost of the HWB example. The computational cost of generating the basis 
functions and AIC is about 250 min in this HWB example. 

CONCLUSION 


A technique for approximating the modal AIC matrix by using basis functions has been 
developed and validated, and a process for using the resulting AIC matrix in aeroelastic analysis 
and design optimization has been proposed. The approximation method has been applied to the 
aeroelastic analyses, and the results are essentially identical to those using direct solution. The 
technique presented has been shown to offer consistent flutter speed prediction on an ATW2 
configuration and a HWB type vehicle with a negligibly small loss in precision. These results 
may have practical significance in the analysis of aircraft aeroelastic calculation and could lead to 
a more efficient design optimization cycle. The basis function approximation approach yields 
significant improvements in computational efficiency as compared to the original approach, 
thereby meeting the objective of this study. 
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TABLES 


Table 1. Summary of total weight and center of gravity location for the ATW2 with different 
fictitious mass configurations. 


Mass configuration 

1 

2 

3 

Target 

Total weight, lb 

4.046 

5.250 

6.450 

2.850 

X center of gravity, in 

14.510 

14.950 

15.220 

13.710 

Y center of gravity, in 

-8.760 

-9.010 

-9.170 

-8.290 

Z center of gravity, in 

0.000 

0.000 

0.000 

0.000 

Number of modes 

15 

15 

15 

10 


Table 2. Summary of natural frequencies (Hz) of the ATW2 with different fictitious mass 
configurations. 


Mode 

Configuration 1 

Configuration 2 

Configuration 3 

Target 

1 

18.10 

15.98 

14.45 

21.32 

2 

62.81 

61.91 

60.60 

63.55 

3 

85.87 

77.10 

71.84 

103.18 

4 

126.60 

112.29 

101.86 

147.79 

5 

191.68 

173.89 

160.47 

208.55 

6 

230.10 

220.25 

215.18 

265.53 

7 

357.83 

317.80 

290.89 

408.75 

8 

403.33 

383.13 

367.52 

434.86 

9 

450.24 

444.21 

426.69 

472.85 

10 

532.04 

477.94 

448.37 

589.94 

11 

628.72 

616.37 

600.44 


12 

633.07 

630.08 

628.51 


13 

708.57 

707.15 

684.02 


14 

729.18 

726.37 

705.74 


15 

782.46 

728.68 

725.76 



Table 3. Summary of the ATW2 flutter results comparison for BFA and direct method. 


BFA method 

Direct method 

Error, % 

Flutter speed, KEAS 

543.59 

543.59 

0.0 

Flutter frequency, Hz 

49.01 

49.01 

0.0 
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Table 4. Computational cost comparison of BFA and direct method. 


Test case 

ATW2 

HWB 

Number of nodes 

269 

-24000 

Number of structural modes 

10 

40 

Number of reduced 
frequencies 

16 

14 

Direct method elapsed time 

43 s 

50 min 

BFA method elapsed time 

8 s 

9 min 


Table 5. Different fuel weight configurations versus 
HWB. 

target fuel weight configurations of the 

Mass configuration 

1 

2 

3 

Target 

Total weight, lb 

403781 

478941 

591681 

535311 

X center of gravity, in 

1337.60 

1334.51 

1331.34 

1332.76 

Y center of gravity, in 

0.75 

0.46 

0.16 

0.29 

Z center of gravity, in 

20.89 

21.35 

21.82 

21.61 

Number of modes 

50 

50 

50 

40 


Table 6. Natural frequencies (Hz) of the HWB with various fuel weight and target fuel weight 
configurations. 

Mode 

Configuration 1 

Configuration 2 

Configuration 3 

Target 

1 

1.60 

1.54 

1.47 

1.50 

2 

1.74 

1.64 

1.53 

1.58 

3 

3.52 

3.35 

3.13 

3.24 

4 

3.64 

3.45 

3.22 

3.33 

5 

5.24 

5.00 

4.54 

4.83 

6 

5.26 

5.04 

4.54 

4.89 

7 

6.42 

5.42 

4.68 

4.92 

8 

6.42 

5.43 

4.76 

4.92 

9 

6.92 

6.66 

5.81 

6.31 

10 

7.12 

6.68 

5.82 

6.31 

11 

7.33 

6.95 

6.11 

6.39 

12 

7.65 

6.95 

6.32 

6.48 

13 

8.22 

7.03 

6.69 

6.85 

14 

8.23 

7.40 

7.06 

7.22 


Table 7. Summary of the HWB flutter results comparison for BFA and direct method. 


BFA method 

Direct method 

Error, % 

Flutter speed, normalized 

0.991 

1.000 

0.86 

Flutter frequency, Hz 

2.671 

2.667 

0.15 
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FIGURES 
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Figure 1. Object-oriented multidisciplinary design, analysis, and optimization tool. 



Figure 2. Flowchart of the flutter analysis module in object-oriented MDAO tool. 
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Figure 3. Aerostructures test wing 2. 
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Figure 5a. Mode 1 (First bending), 21.32 Hz. 
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Figure 5b. Mode 2 (Second bending and first torsion), 63.55 Hz. 

Figure 5. Structural finite element model mode shapes and frequencies of the ATW2 with target 
mass configuration (direct method). 
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Figure 5c. Mode 3 (Second bending and first torsion), 103.18 Hz. 
Figure 5. Concluded. 



Mode 




Mode 


Mode 



Figure 6. Structural finite element model mode shapes error of the ATW2 with target mass 
configuration. 
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Fuel tank area 
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Figure 9. Hybrid wing body aircraft finite element model. 




Flexible mode 1 (First symmetric bending) 


Flexible mode 2 (First anti-symmetric bending) 



Flexible mode 3 (First symmetric torsion) 



Flexible mode 4 (First anti-symmetric torsion) 
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Figure 10. Hybrid wing body mode shapes. 
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Figure 11a. V-g plot. 
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Figure 1 lb. V-co plot. 


Figure 11. V-g and V-co plots for the HWB at Mach 0.50 using direct modal AIC matrices. 
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Figure 12a. V-g plot. 



Figure 12b. V-coplot. 

Figure 12. V-g and V-co plots for the HWB at Mach 0.50 using approximate modal AIC matrices. 
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APPENDIX A 

BASIS FUNCTION APPROXIMATION 

' n \ ( n . N 

Qij= IPkYk A tP J k ¥ k 

K k= ' J V A = I J 

( „ \ 
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V&=1 

f „ \ 

= 'LPkV'k ( A AjVl + A P J 2V2 + --- + A P J nVn) 

\k= 1 

f n \ 
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\k= l y 
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